Terahertz response of dipolar impurities in polar liquids: On anomalous dielectric 

absorption of protein solutions 



ON 

o 
o 

(N 

-t— > 
O 

O 



43 

Oh! 



Dmitry V. Matyushov 

Center for Biological Physics, Arizona State University, PO Box 871604, Tempe, AZ 85287-1604^1 

A theory of radiation absorption by dielectric mixtures is presented. The coarse-grained formu- 
lation is based on the wavevector-dependent correlation functions of molecular dipoles of the host 
polar liquid and a density-density structure factor of the positions of the solutes. A nonlinear de- 
pendence of the absorption coefficient on the solute concentration is predicted and originates from 
the mutual polarization of the liquid surrounding the solutes by the collective field of the solute 
dipoles aligned along the radiation field. The theory is applied to terahertz absorption of hydrated 
saccharides and proteins. While the theory gives an excellent account of the observations for sac- 
charides without additional assumptions and fitting parameters, experimental absorption coefficient 
of protein solutions significantly exceeds theoretical calculations within standard dielectric models 
and shows a peak against the protein concentration. A substantial polarization of protein's hydra- 
tion shell is required to explain the differences between standard theories and experiment. When 
the correlation function of the total dipole moment of the protein with its hydration shell from 
numerical simulations is used in the present analytical model an absorption peak similar to that 
seen is experiment is obtained. The result is sensitive to the specifics of protein-protein interactions 
in solution. Numerical testing of the theory requires the combination of terahertz dielectric and 
small-angle scattering measurements. 

PACS numbers: 77.22.-d, 61.20Gy, 61.25Em, 87.15H-, 87.15np 

Keywords: Dielectric response, mixtures, protein electrostatics, terahertz spectroscopy, structure factor. 
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I. INTRODUCTION 

Dielectric spectroscopy of mixtures is a well- 
established technique which requires theoretical model- 
ing for the data interpretation. The models of dielectric 
response of mixtures traditionally operate by assuming 
that a mixture can be separated into macroscopic di- 
electric bodies. Among the commonly used models are 
the Maxwell- Wagner theory [l[ and various formulations 
of the effective- medium approximation [2]. Both assume 
that a dielectric constant can be assigned to each compo- 
nent, and the latter also requires that the physical prop- 
erties of the host and the impurity are not dramatically 
different. 

The recent rapid development of dielectric techniques 
to study mixtures H, in particular in the terahertz (THz) 
frequency window [J], aims at a different length-scale. 
The interest is mainly driven by the desire to learn about 
electrostatics of nano-scale objects, such as biopolymers 
d IE H H ©, El EH, nano-crystals [nj|, and nano- 
confmed fluids 13]. In particular, one hopes that the 
properties of the nano-scale interface between the sol- 
vent and the solute can be effectively probed by the di- 
electric response. This goal is complicated by the fact 
that essentially any relaxation event linked to electrical 
dipoles in the system contributes to the integral experi- 
mental signal, and theory is required to separate different 
components. While fully atomistic models will be the 
ultimate goal of the theory, it is still useful to develop 
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coarse-grained approaches employing the length-scale in- 
termediate between macroscopic dimensions of classical 
theories [j], Q and fully atomistic length-scale. 

This paper presents a coarse-grained model of the di- 
electric response of dipolar mixtures, aiming in partic- 
ular at the THz frequency window. The model does 
not assume that solutes can be described as dielectric 
bodies, neither docs it assume dielectric continuum for 
a polar solvent. The polar liquid is characterized by 
its wavevector-dependent correlation functions (l4| . and 
a similar approach is invoked for the solutes character- 
ized by their density structure factor. However, instead 
of using completely atomistic structures, the solutes are 
modeled by effective spheres characterized by dipole mo- 
ments, polarizabilities, effective radii, etc. The assump- 
tion of solute sphericity does not pose a fundamental re- 
striction on the theory since it can be extended to so- 
lutes of non-spherical shapes made by overlapping vdW 
spheres of the composing atoms [15[ . However, this sim- 
plification allows us to come up with a set of compact 
analytical equations applicable to analyzing experimen- 
tal data. 

The theory is applied to the analysis of the absorption 
coefficient of THz radiation. Recent measurements on 
hydrated saccharides [l6j] and proteins [l7| have shown 
qualitatively different types of dependencies of THz di- 
electric absorption on concentrations of these two types 
of solutes. The current theory gives an excellent account 
of the observations on saccharides, but fails to repro- 
duce the protein experiments when the dipole moment 
of the protein is assigned to the solute. It is suggested 
that hydrated proteins introduce solvation electrostat- 
ics qualitatively different from the dielectric response of 
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FIG. 1: Schematic diagram of a conventional dielectric 
impedance experiment. The electric field En (in the absence 
of dielectric) is perpendicular to the plane of the liquid film 
such that the field in the dielectric is En/e. This is a longi- 
tudinal field as it sets up the direction of symmetry breaking 
in the homogeneous liquid. The polarization of the cavity 
in the liquid induces the depolarization field and the cavity 
dipole M c opposite to the direction of the external field. The 
average solute dipole (mo), aligned along the external field, 
enhances both the dielectric response and the depolarization 
field of the empty cavity. 



typical dipolar mixtures [la ] . Specifically, hydration lay- 
ers nearest to the protein (ca. 15 A in thickness) be- 
come polarized and thus carry a significant dipole mo- 
ment with the relaxation dynamics different from that of 
the protein This "elastic ferroelectric bag" [18| sur- 
rounding the protein significantly enhances the effective 
dipole moment of the solute observed on the large wave- 
length of THz radiation and can account for the observed 
anomalous dielectric absorption of protein solutions [ItJ ■ 
Since both the scenario of the rigid protein dipole and 
the dipole dressed by the ferroelectric bag can be intro- 
duced into the formalism, the present theory provides a 
tool to separate this new physics from what can be de- 
scribed within the traditional understanding of dipolar 
liquids and solvation electrostatics. 



II. DIELECTRIC RESPONSE OF MIXTURES 

We consider a polar liquid with dipolar impurities (so- 
lutes). The impurities are larger than the molecules of 
the host liquid in most cases of practical interest and 
can physically be realized as molecules or small colloids 
(nanoparticles) . The physics of the problem is clearly 
presented by separating the process of inserting the im- 
purities into two steps: (i) the creation of a hard-core 
cavity in the liquid and (ii) the polarization of the host 
polar liquid by the partial charges of the overall neutral 
solute. We note in passing that the restriction of the 
neutral solute can be lifted when ionic conductivity is 
not an issue, such as the case for many THz dielectric 
measurements. 

The creation of a cavity in a polar liquid results, in 
terms of standard dielectric theories [2(| [2l[ , in a depo- 
larization field, i.e. charges on the cavity's surface that 



create the cavity dipole moment M c opposite to the di- 
rection of the external field. In the standard setup of the 
dielectric spectroscopy experiment shown in Fig. Q] the 
electric field is longitudinal, i.e. parallel to the direction 
of breaking the isotropic symmetry of the liquid by an 
external perturbation. The dipole of a spherical cavity 
of volume fl is then [H| M c = -3P L f2 /(2e + 1), where 
e is the dielectric constant of the homogeneous liquid and 
P L is the longitudinal (superscript "L" ) polarization field 
created by the external source of the electric field Eq. 
Since P L = (e — l)E /(47re), the cavity dipole decreases 
with increasing e. Standard low-frequency (high e) di- 
electric measurements of polar liquids are therefore fairly 
insensitive to impurities. 

The dipole moment of the solute orients itself in the ex- 
ternal field amplifying the dielectric response. This effect 
is partially compensated by an additional polarization of 
the cavity surface by the internal dipole acting to enhance 
the cavity dipole in the direction opposite to the exter- 
nal field (Fig. [J) . The solute dipoles can be considered 
as independent in the limit of infinite dilution, and the 
change of the dielectric response is linear in the dipoles' 
concentration. This approximation limits the range of 
concentrations by the requirement that the Onsager ra- 
dius of the solute-solute dipolar interactions is below the 
average distance between them. 

The situation becomes more complex for a finite con- 
centration of solute dipoles when an additional effect of 
their collective field gains in importance. The alignment 
of solute dipoles in the external field creates a net average 
dipole moment (mo) (Fig. [T]) and a corresponding non- 
zero net electric field that can potentially polarize cav- 
ities and alter their cavity dipoles. Since internal fields 
are commonly large compared to the external field, this 
effect, nonlinear in the solutes' concentration, can be po- 
tentially significant. 

The arguments we have presented so far apply to the 
standard dielectric impedance measurements employing 
longitudinal electric fields. THz experiments employ a 
different geometry where the absorption of a pulse of 
electromagnetic wave propagating orthogonally to a thin 
(ca. 100 /mi) film is measured [6, 9]. In this case, the 
electric field is transversal, i.e. it is perpendicular to 
the direction of axial symmetry breaking introduced in 
the isotropic liquid by the direction of the wavevector 
[22l I23I ]. One measures then the transverse dielectric 
response and the cavity is polarized differently. The 
dipole moment of the cavity along the field becomes 
M^ 1 = — P T Oo x 3e/(2e + 1), where the transverse po- 
larization (superscript "T") is P T = (e - l)E /(47r). It 
is clear that the cavity dipole produced in response to 
the transversal field is not screened by the high dielectric 
constant of a polar liquid. Microwave absorption mea- 
surements are therefore expected to be significantly more 
sensitive to impurities than conventional dielectric mea- 
surements. This distinction is the physical basis of the 
sensitivity of the transversal absorption experiments to 
electrostatic changes in molecular or nano-scale solutes 



3 



in the following form 



III. RESPONSE FUNCTION 

We now turn our attention to a detailed analysis of the 
transverse dielectric response of dipolar mixtures. In or- 
der to approach this problem we will use the approxima- 
tion of linear response of the solvent to the electric held of 
the solute. The linear response approximation states that 
the solvent response function is insensitive to the magni- 
tude of the solute electric field and in fact can be calcu- 
lated for a fictitious solute with all partial charges turned 
off (zero dipole for a dipolar solute) [l5[ ■ Even though the 
electrostatic response is linear, the response to the solute 
repulsive core cannot be calculated within linear mod- 
els since the repulsive potential of the solute produces 
a large and nonlinear perturbation of the solvent struc- 
ture. This perturbation renormalizes the spectrum of 
the solvent fluctuations modifying the linear (Gaussian) 
response function (24|. In dielectric theories, this modi- 
fication is included by imposing boundary conditions on 
the solution of the Poisson equation. The problem be- 
comes way more complex at the molecular level and is 
commonly solved in terms of angular-dependent distri- 
bution functions (25J. 

We will adopt here Chandler's formulation of the Gaus- 
sian model [2J] in which the linear response function, 
modified by the presence of solute, is sought by impos- 
ing the condition of vanishing solvent density from the 
solute's hard core. In case of polarization response, this 
condition implies the polarization field P vanishing from 
the hard core of the solute. One can then define a gen- 
erating functional of the polarization field as follows [l4[ 

G[E ] = J exp [-(/3/2)P * x^ 1 * P + /?E * P] 



Jp[P(r)]OT. 



(1) 



Here, Eo is an external electric field, the asterisk denotes 
both the volume integration and tensor contraction, and 
(3 = l/(k^,T) is the inverse temperature. Further, \ s 
is the 2-rank tensor of the Gaussian fluctuations of the 
polarization field in the homogeneous solvent and the 
product of delta functions runs over all points within so- 
lute's hard-core of volume f^o and over all solutes (index 
i) . This term ensures that the polarization field vanishes 
from the volume of each solute in the mixture. 

Functional derivatives of G[E*o] over the external field 
Eo produce correlation functions of the polarization field 
of the solvent in the presence of solutes. The Gaussian 
integral over the polarization field P(r) can be calculated 
exactly resulting in a Gaussian functional in the external 
field Eo. The corresponding renormalized response func- 
tion x gains most compact representation in the inverted 
k-space [l4j . It can be written in the k, u>- representation 



x(ki,k 2 ,u;) = Xs(ki,w)5 kl! k 2 

- ]T X R (ki,u;) ■ e*( kl - k2 )- r %(ki - k 2 ) • Xs (k 2) u). 

i 

(2) 

Here, <5ki,k 2 = (27r) 3 (5(ki — k 2 ) and #o(k) is the Fourier 
transform of the step function defining the excluded vol- 
ume of the solute the translational dynamics of which 
are neglected. The direct-space Heaviside function #o( r ) 
is equal to unity within the solute and is equal to zero 
outside the solute. The inverted-space function is given 
by the Fourier transform 



?o(k) 



(3) 



where integration is over the solute volume fio- 

The response function of the mixture x(ki,k 2 ,w) de- 
pends on two wavevectors ki and k 2 separately, instead 
of ki — k 2 of the homogeneous liquid, because of the in- 
homogeneous response produced by each solute marked 
by index i. This response function combines the dipo- 
lar response function of the homogeneous liquid X s (k, lu), 
the information about the solute shape incorporated into 
#o(k), and the renormalized function x^k, u>) (see be- 
low). 

The response function of an axially-symmetric dipolar 
liquid is expandable into longitudinal (L) and transverse 
(T) projections [H, HI 



(4) 



where J L = kk and J T = 1 — kk are the orthogonal 
longitudinal and transverse dyads. The k = values 
of the response projections are directly related to the 
frequency-dependent dielectric constant of the host liquid 



47rxf(0,w) 
4^(0, w) 



e(o,) - 1. 



(5) 



The entire k, w-dependence of the projections x L ' T {k 1 uj) 
is given in Ref. 1271 but only the transverse projection is 
required for the problem considered here (see below). 

The last function in Eq. ^ that requires definition is 
X R (k, L)). This function appears in the solution for the 
generating functional in Eq. ([T]) as a result of renormaliz- 
ing the dipolar response of the homogeneous liquid by the 
solute cavity. It thus contains the information about both 
the solvent and the solute [13, H3] ■ Only k = transverse 
projection of this function appears in the equations for 
the transverse dielectric response of the dipolar mixture 
and that is given by the following equation 



X R ' T (0,u,) 



3e(cj) 
2e(w) + 1' 



(6) 



We will now use Eq. ^ to calculate the transverse 
dipole moment M t {lo) of the dielectric sample produced 
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in response to the electric field of the electromagnetic 
radiation oscillating with frequency ui 



V (t)=e T E e l 



(J) 



Here, the polarization unit vector e T is perpendicular to 
the direction of propagation k. 

The dipole moment M T (u>) combines two contribu- 
tions: the dipole moment induced directly by the exter- 
nal field of the radiation (radiation wave-length is much 
larger than any molecular scales in the system) and an 
additional collective polarization induced by the solute 
dipoles aligned along the external field. These two con- 
tributions are described by correspondingly the first and 
the second summands in the following relation 

M t (lo) = e T ■ x(0, 0, u) ■ e T E 

+ e T ■ X (0, T ( k K k - ri ' mo», (8) 

i 

where T(k) is the Fourier transform of the dipolar tensor 
T = — V r V r '|r — r'| _1 and, as above, the asterisk refers 
to integration over k-space and tensor contraction. In 
addition, moi(o>) is the solute's dipole moment aligning 
itself with the oscillating external field. 

The dipole moment moi(w) is a sum of two compo- 
nents: the electronic dipole induced instantaneously (on 
the time-scales of interest) by the external filed and a per- 
manent dipole inertially rotated by the torque imposed 
by the external field. The inertial component can be cal- 
culated from the linear-response approximation [25[ with 
the result 

e T • m 0l (w) = (a 0e + aj£„[l - iu>$(-w)]) f d {u)E - (9) 

Here, ao ie is the solute electronic polarizability and the 
permanent dipole polarizability is given as 



(10) 



where is the transverse Kirkwood factor [28| of the 
correlated orientations of the solute dipoles: 



T 

9a. K 




ej - (§j ■k)(k-e 3 ) 



(11) 



Here, ej are the unit vectors of the solute dipoles. If the 
dielectric constant eo can be assigned to the solutes, then 
e = 1 + 2Trf3mlp gQ K , po = N /V. We also assumed 
isotropic polarizability of the solute and, in addition, for 
solution problems, the permanent dipole too should be 
properly renormalized from the g as-p hase value by the 
effect of the solute polarizability |2ll . l29l [30| . Further, 
the factor fd(u >) in Eq. (|5J) is the Onsager directing field 
correction [2 lL |29| accounting for the difference between 
the electric field of the radiation and the local electric 
field imposing torque on the solute dipole. 



The Laplace-Fourier transform $(w) in Eq. ([9]) repre- 
sents correlated rotational dynamics of the solute dipole 



e 1 • (mo(t)Mo(O)) • e T e iut dt 



(12) 

where M = J2j m o,j 1S the total solute dipole in the 
sample. In case of a single-time Debye rotational relax- 
ation with the relaxation time To the term in the square 
brackets in Eq. © gains the form 



1 - iw$(—w) = (1 + iwT )" 



(13) 



The Debye approximation in Eq. (|13p is typically suf- 
ficient for rigid dipoles dissolved in a polar solvent. 
The situation potentially becomes more complex for soft 
nano-scale solutes, biopolymers in the first place. The 
dynamics of the dipole moment is then affected by low- 
frequency vibrations 0] altering id). As we discuss 
below, the inclusion of a non-vanishing dipole moment 
of the protein's hydration shell, with its own dynamics, 
makes the problem even more non-trivial, further com- 
plicating the form of $(— u>). 

The first term in Eq. ([8]) can be easily calculated by 
combining Eqs. ([5]) and © and noting that 6 (0) = £1 
[Eq. This calculation then results in a simple rela- 
tion for the difference between the response function of 
the mixture Xmix(^) = M t (uj)/(VEo) and of the homo- 
geneous liquid 47rx(w) = e(u>) — 1 



4?rAx = -r? /(w). 



(14) 



Here, Ax = Xmix(^) - x( w X = N n flo/V is the volume 
fraction of the solutes in the mixture with the overall 
volume V, and 



3eM(eM - 1) 
2e(w) + 1 



(15) 



We note here that a more simple (and elegant) deriva- 
tion of the response function of a low-concentration mix- 
ture as given by Eqs. (Tj"4"l) and ([15)) can be found in 
Ref. Equations (Ti"4")) and (fl5]) also represent a low- 
concentration limit of the Maxwell- Wagner formula [l[. 
Our microscopic consideration is thus consistent with 
macroscopic arguments. The microscopic description is 
however required to get correctly the second summand in 
Eq. ([5]) describing the collective response of an ensemble 
of solute dipoles. This is what we consider next. 

The response function of the solvent to the presence 
of the solute includes two parts corresponding to two 
summands in Eq. The first summand represents 

the response of the liquid to an infinitely small solute 
which does not perturb the spectrum of dipolar fluctua- 
tions of the liquid. This part is easy to calculate and its 
relative contribution to the response is A^(w)/x(w) = 
y (oj)f d (uj), where 

2/oM = (47r/3)p ao,e + (27r/3)/3g r K TOoPo [1 - iw$(-cj)] 

(16) 
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is the dipolar density of solutes defined in analogy with 
a similar quantity of homogeneous liquids [l| . 

The contribution from the second term in Eq. {2} is 
the correction of the solvent response introduced by the 
excluded volume of the solute. This calculation is more 
complex. After some algebra one arrives at the mixture 
susceptibility relative to the susceptibility of the homo- 
geneous liquid 



-Vo; 



3e(w) 



2e(w) + 1 
+ 2/oM/dM I 1 



3e(w) 
2c(uj) + 1 



I(w,r)o,R) 



(17) 



The only non-trivial part in this equation is the inte- 
gral I(aj,r]0i R), arising from the combined effect of the 
volume excluded by the solute from the solvent, many- 
body solute-solute correlations, and microscopic correla- 
tions between the dipoles of the solvent. It is given by 
the relation 



I(u,Vo,R) 



QR e(0) - 1 
7r e(to) — 1 



dkjf(kR)S a (k,r) ) 



(18) 

in which j% (x) is the first-order spherical Bessel function 
and R = (ero +cr) /2 is the distance of the closest approach 
of the water molecules with the effective hard-sphere di- 
ameter a to the solute characterized by its hard-sphere 
diameter do- 

The density-density structure factor So(k,rjo) in Eq. 
(fT5)) is responsible for a nonlinear dependence of the re- 
sponse function of the mixture on the solute concentra- 
tion. The k = value of the structure factor S o (0, 
(So(0, 770) — ► 1 at ?7o — * 0) is the reduced compressibil- 
ity of the solute component of the mixture. It is equal 
to the experimentally measurable osmotic compressibil- 
ityPil 



So{0,Vo) = Xc 



= ( dpo \ 



(19) 



where II is the osmotic pressure and the derivative is 
taken under the condition of osmotic equilibrium. 

The transverse dipolar correlation function xj(fc,aj) 
in Eq. (|T5|) does not depend on the solute concentration, 
but incorporates spacial transverse correlations between 
dipoles in the polar liquid. We provide its functional 
form here for completeness and refer the reader to Refs. 
I27I and [H for a more detailed account of this problem 



1 



xT(k,u) S T (0) l+p'(fca) 2 



6(0) - 1 



1 



(20) 



In this equation, S T (k) is the static structure factor of 
transverse dipolar fluctuations. A simple extension of 
the mean-spherical solution for dipolar fluids [2611 gives 
S T (k) consistent with numerical simulations [33j. This 
formalism is used here for numerical calculations of the 




FIG. 2: I(u>,rio,R) calculated at to = and ui — 1 THz as 
indicated in the plot vs the reduced distance of closest solvent 
approach to the solute R/a = (00/0"+ l)/2. The dashed line 
is the dielectric-continuum result of Eq. (|21|l ; 770 =0.1 and 
the solvent parameters are those of water (see Appendix). 



function I(u>,Vo) in Eq. (|18p. Finally, the parameter p' 
in Eq. ([20]) quantifies the relative contribution of trans- 
lational vs rotational motions of liquid's dipoles in the 
overall response as discussed in Ref. [3J. 

The approximation of continuous dielectric corre- 
sponds to the neglect of the k-dependence in the trans- 
verse response function xj(^w) in Eq. (fl~8f assuming 
xJ(k,Lu) ~ x^(0i w )- The dependence of frequency then 
disappears from the integral I(to,r]o) which simplifies to 



/(TO, R) = (6RM / dk 3 l{kR)S {k, 770). 



(21) 



The dielectric-continuum integral I(j]q,R) is equal to 
unity for an ideal solution when S(k, 0) = 1. This ideal- 
solution/continuum limit then results in a simple equa- 
tion for the mixture's dielectric response 



-Vo 



3e(w) 
2e(w) + 1 



yo{u)fd(u, 



e(w) - 1 



2e(cj) + 1 

(22) 

It shows that the presence of very dilute solute dipoles 
lowers the transverse response because an enhanced de- 
polarization of the cavity wins over the direct alignment 
of the solute dipoles by the external field. It is clear 
that this result cannot sustain itself as the concentration 
of dipolar impurities grows since the limit of a negative 
dielectric constant can potentially be reached. Solution 
non-ideality must slow the negative decay of the mixture 
susceptibility or change its sign to positive. 

The continuum integral I(r]o,R) can be rewritten in 
r-space as 

1(770, R) = 1 + (Po/th) / dridr 2 /o(ri)/io(ri 2 )/o(r2), 



(23) 

where /io( r i2), ^12 = fa — r i is the pair correlation func- 
tion of the solutes and /o(r) are Mayer /-functions rep- 
resenting hard cores of the solutes. If long-ranged in- 
teractions between the solutes are neglected, the lowest- 
order density expansion of the pair correlation function 
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^-o( r i2) yields the third virial coefficient C\n of the mix- 
ture of hard spheres of diameter R (components 1) and 
diameter oq (component 2): 1 = 1 — (3po/^o)Cn2, 
fio = (4-7r/3)i? 3 . The third virial coefficient of the hard- 
sphere mixture is known [35l . l36l ] . For solutes much larger 
than solvent, one can put R ~ 00/ 2 with the result 
1(770) = 1 - (W 8 )( 3 + 65/24). This simple equation 
compares reasonably well with the direct numerical inte- 
gration using the Percus-Yevick (PY) density structure 
factor. The numerical integrals can be approximated by 
a polynomial of R/<Tq and r] , and this fit is provided in 
Appendix for < 770 < 0.3. 

The interactions between hydrated proteins are com- 
plex, and the structure factor from hard-sphere repul- 
sions can be used only for a limited range of ionic 
strengths when Coulomb forces are sufficiently screened 
[37l [3o | . The structure factor Sn(k, tjq) is directly mea- 
sured by small-angle scattering [38l l39j | and can be nu- 
merically reconstructed from a linear combination of a 
repulsive and attractive potentials; a combination of 
Yukawa potentials is often used [3{|. The small-fc part 
of the structure factor is strongly affected by long-range 
interactions, and there is a peak at q m ~ 27r^/no at the 
average distance between the solutes in solution. Since 
the amplitude of the peak is relatively small [H[ , a gen- 
eral insight into how correlations between hydrated pro- 
teins affect the dielectric response can be gained from 
an empirical approximation for So(k,rjQ). The following 
approximation (analogous to the empty core model (40| ) 
follows directly from the low-density expansion of the di- 
rect correlation function of a hard sphere 

5 (fc,77o) = [l + oji(fc«7o)/(*£ro)] -1 (24) 

in which the constant a is chosen to reproduce the os- 
motic compressibility a = 3(<So(0, rjo)" 1 — 1). The re- 
sulting integral is just a function of a (R ~ 0o/2). Its 
numerical value can be approximated by a Pade form, 
I{a) = (l + 0.0908308a-0.00226567a 2 )/(l + 0.131266a- 
0.00434023a 2 ), which allows one to use the osmotic com- 
pressibility, affected by both repulsions and attractions, 
as input to obtain the dielectric response. The approx- 
imation in Eq. (|24p is accurate up to 770 ~ 0.1 when 
compared to the direct integration with the PY density 
structure factor. 

It is worth noting at this point that the continuum ap- 
proximation is inaccurate at low frequencies u> ~ over- 
estimating the cavity polarization in the entire range of 
solute sizes of common interest (Fig. [5]). This happens be- 
cause of a very sharp decay of the structure factor S T (k) 
at small k- values [33j which, in the continuum limit, is 
replaced by its k — value S T (0). The continuum ap- 
proximation becomes more accurate as the frequency in- 
creases and the dielectric constant drops (Fig. but it 
needs to be tested before applied in a frequency range of 
interest. Nevertheless, in the range of THz frequencies, 
the continuum limit [Eq. Q2ip] presents a useful simplifi- 
cation of Eq. (TTS|) . which, in conjunction with Eq. (pM|) . 
yields the dielectric response solely in terms of observable 
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FIG. 3: Relative change in the absorption coefficient of the 
mixture at uj — [see Eqs. (|26[) and (|27|l ] as a function of 
the volume fraction of the solute 770. The solid line refers to 
the entire response function from Eq. (|17[) , while the dashed 
line shows the contribution of the first term only. The latter, 
linear in 770, is the limit of zero solute dipole. The overall non- 
linear dependence on 770 is the result of mutual polarization of 
the cavities by the solute dipoles. The solute size and dipole 
are those of the A6-85 protein. 



quantities. 

The actual dependence of the dielectric response on 
the solute volume fraction is more complex than a nearly 
linear decay suggested by Eq. (f2"2"]l . It is shown in Fig. [3] 
where a static w = response is calculated for parame- 
ters specific to A5_86 protein discussed below. Baxter's 
solution of the PY closure [25[ for So(k,r] ) was used in 
these calculations. In addition, the microscopic trans- 
verse response function of the solvent dipoles was taken 
according to Ref. [2?] and the static structure factor was 
calculated from a corrected mean-spherical approxima- 
tion suggested in Ref. HI. The dependence of the dielec- 
tric response on 770 is curved down, thus eliminating the 
dielectric catastrophe following from the linear extrapola- 
tion of Eq. (|2"2")) . However, the shape of the concentration 
dependence depends on frequency, and the curvature is 
just the opposite one for the THz response (see below). 

A notion regarding theory's approximations is relevant 
here. One might argue that the point-dipole model is 
too restrictive for the electrostatic field of a protein with 
typically a non-zero overall charge and the prevalence of 
charged residues on its surface. We believe that the ap- 
proximations adopted here are adequate, and the theory 
might actually be more quantitative than it seems. First, 
the solvent response function is independent of the solute 
charge in the linear response approximation [33| and is 
identical to the one obtained for a fictitious solute with 
all charges turned off. The linear response approximation 
might obviously fail, and that certainly puts a restriction 
on the current theory. Second, the perturbation Hamil- 
tonian for the current problem is the interaction of the 
sample dipole moment with the external electric field of 
the radiation. Since the THz wavelength obviously ex- 
ceeds any molecular dimension, a dipolar approximation 
is appropriate for solutes of nano-scale dimension. Fi- 
nally, the total solute charge can contribute to conduc- 
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tivity [4l[ that is normally subtracted from the dielectric 
response and is insignificant in the THz frequency range. 
The dipole moment of charged solutes is then defined 
relative to the solute's center of mass 1411 . 



IV. COMPARISON TO EXPERIMENT 

One of parameters reported in THz dielectric measure- 
ments is the relative absorption coefficient Aa(w)/a(u), 
where Aa(u>) = a mix (uj) — a(uj) is the change in absorp- 
tion coefficient of the mixture relative to the pure liquid. 
The absorption coefficient is defined [42|, as the ratio 
of the rate of energy dissipation by the medium (£}u over 
the Poynting vector S(lu) of the incident radiation 



a(ui) 



(£)» 



(25) 



By combining the standard equations for the Poynting 
vector in dielectric media [20L l43j with energy dissipation 
in terms of the dielectric response function one gets 
the equation 



a(u>) 



4-7TW x"M 



(26) 



which can be applied either to the mixture or to the pure 
liquid (c is the speed of light in vacuum) . 

Assuming that the deviation of the response Ax(u>) 
caused by impurities is small compared to the dielectric 
response of the pure liquid, one can easily derive an ex- 
pression for the relative change of the absorption coeffi- 
cient 



Aa(w) 4ttAx"M 2ttAx'(w) 



a(u>) 



(27) 



In this equation, the variation of both the imaginary and 
the real parts of the response are taken into account when 
impurities are introduced into the polar liquid. In partic- 
ular, for solutes with small dipole moment, one can drop 
the term proportional to yo(to) in Eq. (|17p and arrive at 
a simple relation 



Aa(w) 
a{w) 



/» /» 



'H 2e'(w) 



(28) 



where f(u>) is given by Eq. (flU 

Figure |4] shows the comparison of Eq. (|28|) (lines) with 
the experimental dependence (points) of the absorption 
coefficient on the concentration of trehalose dissolved in 
liquid water [T^]. The details of the calculations and 
the parameters used to produce the plot are given in the 
Appendix. Because of the small dipole moment of tre- 
halose, a complete calculation of the dielectric response 
function of the mixture is not required (the term pro- 
portional to yo (ui) in Eq. (fTT)) is small) and Eq. (|28p is 
sufficient. The dashed and dash-dotted lines in Fig. Q] 




FIG. 4: Aa(w)/a(u) at u = 2.5 THz calculated from Eq. (|28|) 
(solid line) and measured experimentally [16| (points) for the 
aqueous solution of trehalose. The dashed and dash-dotted 
lines represent, correspondingly, contributions from the first 
and second terms in Eq. (|28[) such that the solid line is their 
difference. The frequency-dependent dielectric constant of 
water in the THz range of frequencies was taken from Ref. 
and the molecular volume of trehalose f2o = 278 A 3 [45l ] 
was used to convert from experimentally reported molar con- 
centrations to volume fractions (see Appendix for the details 
of calculations) . 




0.02 0.03 
^0 



0.05 



FIG. 5: Relative change in the absorption coefficient at 2.25 
THz mimicking the solution of protein A6-85 studied by THz 
spectroscopy in Ref. Il7l . The solid line is the dielectric re- 
sponse calculated from Eq. (|17[) with mo = 61 D ((/Jk = 2/3) 
and the effective radius of 12.1 A and the points are experi- 
mental measurements [l7| . The dashed line refers to the first 
term in Eq. (|17|l representing the polarization of the solute 
cavities by the external electric field. 



show the first (imaginary part) and second (real part) 
terms in Eq. (f2"8")) . It is clear that changes in the imagi- 
nary and real parts of the dielectric susceptibility upon 
the addition of impurities are comparable in magnitude 
and should both be included. The only solute parame- 
ter entering Eq. (|28p is its volume. Equation Eq. |28p 
can therefore be used to determine molecular volumes 
of weakly polar solutes by means of dielectric measure- 
ments. 

In an attempt to see what might be the theory predic- 
tion for the case of protein solutions we have mimicked 
the conditions reported in Ref. [13 where the absorbance 
of the solution of a five helix bundle protein A 6 _ 85 (4r| 
showed a maximum at the volume fraction of protein 
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FIG. 6: Change in absorption coefficient of solution of plasto- 
cyanin relative to bulk water at 2.25 THz. Dipole correlation 
function of the protein with the surrounding water shell was 
taken from MD simulations The solute volume fraction 
7/0 was calculated by adding the width of the hydration shell 
(20 A) to the radius of plastocyanin (16.8 A). The points 
represent the experimental results for A6-85 protein [TtJ re- 
calculated from experimental molar concentrations by using 
the combined volume of plastocyanin and its polarized water 
shell. The solid curve refers to the calculations done with the 
hard-sphere solute-solute structure factor in Eq. (fT8|) . In or- 
der to show the sensitivity of the results to the solute-solute 
correlations, the dashed line represents the continuum inte- 
gral 7(r/o) = 1 — 3t/o with the slope against 770 much exceeding 
that for hard-sphere solutes. The inset shows the change in 
dielectric loss of the plastocyanin solution at 770 = 0.05 rela- 
tive to bulk water against frequency measured in 10 12 s _1 . 



below 1% (points in Fig. [5J. The calculations (see Ap- 
pendix for the parameters used) show almost no effect 
of proteins' dipoles and a negative contribution to the 
absorption, as in the case of trehalose above and in an 
obvious disagreement with the experiment. 

There is also a clear difference between Figs. [3] and [5] 
While Fig. [3] shows a clear effect of the mutual cavity 
polarization by solutes' dipolar fields for the same set of 
parameters, there is almost no effect of the solute dipolar 
component in Fig. [5] (cf. solid and dashed lines). The 
difference comes from the dynamical effect. The solute 
dipoles do not have time to reorient on the time-scale 
of the THz pulse and the corresponding contribution is 
strongly diminished by the relaxation 1/(u>tq) term. The 
THz pulse thus probes almost exclusively the electronic 
polarizability of the solvated proteins. 

For the solute dipoles to be seen in the THz response, 
either a much faster relaxation or a significantly larger 
effective dipole are required. Faster relaxation of pro- 
tein's dipole seems improbable given that numerical sim- 
ulations show an almost exclusively single-component ro- 
tational relaxation with the relaxation time in the range 
3-6 ns [4l[ • The hydration shell thus emerges as the most 
probable candidate to explain the differences between the 
theory and experiment. 

In order to obtain more quantitative insights into the 
problem, results of numerical simulations of protein so- 



lutions are required. We found recently [H, [l9| that, in 
accord with the suggested interpretation of experimental 
THz data [l7|, proteins are capable of polarizing their 
hydration shells ~ 10 — 15 A into bulk water. This po- 
larization results in a significant non-zero average dipole 
moment of the hydration shell (Im^l), which reached the 
value of ~ 10 3 D in simulations of metalloprotein plas- 
tocyanin [l8j . The dynamics of this ferroelectric cluster 
around the protein are however decoupled from a much 
slower tumbling of the protein occurring on the time-scale 
of nanoseconds. The relaxation of the shell's dipole 
is clearly two-component, with a very fast initial relax- 
ation on a sub-picosecond time-scale, followed by a low- 
amplitude tail lasting hundreds of picoseconds. The fast 
component correlates with low-frequency vibrations of 
the protein deforming water's "elastic ferroelectric bag" 

In this picture, the solute dipole mo should be re- 
placed with the sum of protein's and shell's dipoles 
M = mo + m^, . The dynamics of this total dipole gives 
input to determine function which, together with 

(M 2 ), yields yo(ui) [Eq. (fT?j)) ]. These parameters were 
extracted from simulations of plastocyanin carrying the 
negative charge of —8 in its oxidized state and hydrated 
by N w = 21076 TIP3P waters (THH. The shell of wa- 
ter molecules of width 20 A was added to the effective 
radius of the protein to obtain the effective radius of the 
protein/water cluster and the volume fraction of coupled 
protein/water dipoles in solution (see Appendix for de- 
tails). The dielectric response of the solution was then 
calculated from Eqs. (T7]) and (j2"6]) . 

Figure [6] shows the concentration dependence of the 
solution absorption coefficient with the PY hard-sphere 
structure factor So(k,rjo) (solid line). The points, shown 
for reference, are data on Ag-85 protein (l7| rescaled with 
the volume of the plastocyanin/water cluster. The cal- 
culation indeed yields a maximum in the absorption co- 
efficient which turns to negative values with increasing 
volume fraction. The outcome of these calculations is 
sensitive to the form of the density structure factor and, 
therefore, to protein-protein interactions in solution. In 
order to illustrate this point, the dashed line in Fig. [5] 
shows the result of calculations with a stronger effect of 
repulsions and thus a steeper decay of 6*0(0, 770) with in- 
creasing 770 . 

The hard-sphere model might not be adequate for all 
proteins and electrolytes. For instance, for the ionic 
strength employed in Ref. [l7| (0.05 M), the interactions 
between hydrated bovine serum albumin (BSA) proteins 
are dominated by electrostatic repulsions [47j . These pro- 
teins are negatively charged, similarly to plastocyanin, 
and the long-range interactions are dominated by the 
screened Coulomb potential. The osmotic compressibil- 
ity 5(0, 770) of BSA quickly drops with increasing pro- 
tein concentration to the level 5(0,770) ~ 0.1 — 0.2 and 
then does not significantly change when the concentra- 
tion is further increased [47| • With such a dependence of 
5(0, 770) on the volume fraction 770 the peak in absorption 
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FIG. 7: Change in absorption coefficient of solution of A6-85 
protein at 2.25 THz relative to bulk water. The normalized 
dipole correlation function of the protein with the surround- 
ing water shell 3>(w) was taken from MD simulations of plas- 
tocyanin metalloprotein [lU . The density structure factor 
of the proteins was constructed from the hard-sphere and 
Yukawa effective potential [Eq. (|29[) ] obtained in Ref. [Hj] 
by fitting small-angle scattering data. The curves refer to 
different dipole moments (in D) of the protein-water cluster 
as indicated in the plot. The cluster dipoles were assumed 
uncorrelated, = 2/3. The points are the experimental 
results from Ref. [l7| converted to volume fraction with the 
hard-sphere diameter oo = 37.8 A from the effective protein- 
protein interaction potential [iH ]. 

vanishes (Fig. [5]). Note that no absorption peak against 
protein concentration was detected for BSA in dielectric 
terahertz measurements at ui = 1.56 THz [§[. 

The inset in Fig. [5] shows the frequency dependence 
of the dielectric loss Ae"(u>). As is seen, the change of 
the loss relative to bulk water can be either positive or 
negative, depending on the frequency range. A complex 
concentration dependence seen for the absorption coeffi- 
cient in Fig.[6]is the cumulative effect of the concentration 
dependencies of e' mi ^(ui) and e^jx(w). 

The A6-85 protein is uncharged and the corresponding 
protein-protein interaction can be modelled either as a 
sum of soft repulsion and exponentially decaying attrac- 
tion or, alternatively, as a sum of hard-sphere (mhs) and 
attractive Yukawa potentials [48[ : 

u{r) = u HS (r) - e (a /r) e~^-^' 5 e{r - a ). (29) 

In Fig. [7J we used this latter approximation for the in- 
teraction potential to calculate So(k,r]o) 49] and then 
applied this structure factor to the calculation of the 
THz absorption coefficient. In the absence of dipole mo- 
ment dynamics for this protein, we used the normalized 
self-correlation function of the protein- water dipole from 
plastocyanin simulations A set of curves in Fig. 

[7J refer to different values of the dipole moment of the 
protein-water cluster, with the lowest curve correspond- 
ing to the protein dipole alone. Qualitatively, the absorp- 
tion curves do go through maxima with increasing dipole 
of the solute, and the protein solution absorbs stronger 
than bulk water. However, the maxima are broader than 



in experiment and the agreement is only qualitative at 
best. 



V. DISCUSSION 

The present model of the dielectric response targets 
physical situations when large solutes dissolved in po- 
lar solvents do not extend to dimensions of a dielectric 
material. Large cavities in polar liquids carry depolariza- 
tion dipoles oriented oppositely to the external field, with 
their magnitudes scaling linearly with the solute volume. 
These depolarization dipoles accumulate a negative con- 
tribution to the absorption coefficient. The intrinsic so- 
lute dipoles, which align along the external electric field, 
increase the absorption and also produce a non-zero local 
electric field that re-polarizes neighboring cavities. This 
collective effect, non-linear in the solute concentration, 
is sensitive to the solute-solute correlations and is de- 
scribed by convoluting the solvent dipolar response with 
the density structure factor of the dissolved solutes. 

This model performs exceptionally well when tested 
against experimental THz measurements for weakly po- 
lar impurities (Fig. [5]). In this case, only depolarization 
of cavities contributes to the response, and that part of 
the problem seems to be well captured by dielectric the- 
ories. Even though solvation of saccharides distorts the 
structure of water on the microscopic scale 15 0L 1 5 ill and 
slows down the dynamics of the hydration layer [52| , THz 
absorption seems to be insensitive to such changes, and 
the resulting signal is well described by a purely dielectric 
response. This conclusion is consistent with the recent 
light scattering spectra of trehalose solutions [H[ sug- 
gesting only a local perturbation of the water structure 
restricted to the first solvation shell, which is typical for 
many small molecular solutes. 

Polar impurities introduce both the effect of individ- 
ual solute dipoles and their collective polarization effect. 
The response-function formalism employed here does not 
involve any large-scale changes in the solvent structure 
induced by the solute. This formulation then fails to re- 
produce the anomalous increase in the ab sorp tion of pro- 
tein solutions over that of bulk water 0, [U . Computer 
simulations [l8[ show instead a high extent of coopera- 
tivity between hydration shells and protein's motions. In 
addition, a significant polarization of the water shell ex- 
tending 10-20 A from the protein surface into the bulk is 
observed. When the magnitude and correlation function 
of the protein- water total dipole are substituted into the 
equations for the solution response, the theory shows a 
maximum in the absorption coefficient qualitatively sim- 
ilar to experimental observations. The maximum can 
therefore be considered as an observable signature of the 
"elastic ferroelectric bag" found by simulations [18j . The 
shape of this anomalous absorption maximum is however 
sensitive to the interprotein interaction potential and will 
be affected by several factors including protein's ioniza- 
tion state and the ionic strength of the solution. 
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APPENDIX A: DETAILS OF CALCULATIONS 

The dependence of the absorption coefficient on 
frequency arises predominantly from the frequency- 
dependent dielectric constant of the solvent. Dielectric 
measurements of water [44j in the THz range, extended 
to more typical low-frequency dielectric values, have been 
used to produce Figs. [2] and 0]- [71 The dielectric constant 
is given by the following relation [44| 



TABLE I: Solute parameters used in the calculations. 



Solute 


(ao/2)/ A 


m /D 


ro/ns 


Trehalose 


8.2 


1.75 


0.05 


^6-85 


12.1" 


61 b 


3 


Plastocyanin c 


16.8 d 


248 e 


2.8 / 



"From Ref. |17| . The following set of parameters from Ref. [48ll 
was used to represent the protein-protein interaction potential in 
Eq. J29}: <j = 31.8 A, e/k B = 419 K, S = 4.14 A. 

''Calculated from equilibrated protein geometry and atomic par- 
tial charges l53l . 

c According to MD simulation data from Ref. f 1SH - 

d From the vdW volume of the protein using the Amber FF03 force 
field. 

e (mo) calculated from the MD trajectory relative the center of 
mass, total charge of the Ox state of the protein is —8. Fluctuations 
of the protein dipole are caused by protein's vibrations. 

■'Calculated from the exponential fit of the time self-correlation 
function of the protein dipole. 



Aer 



Ae 2 



three-exponent fit of the simulated correlation function 



1 — lOJTi 1 — ILOT2 



(Al) 

where Ae 1 = 73.9, Ae 2 = 1.56, and £oo = 2.34. The De- 
bye relaxation times and the parameters of the resonant 
component are: n = 8.76 ps, t 2 = 0.224 ps, uj s /2it = 5.3 
THz, 7s /27T = 5.30 THz, A s /{2tt) 2 = 35.1 THz 2 . The 
parameter fd(oS) [Eq. |9])] accounts for the difference be- 
tween the external and the local directing (torque) fields. 
It depends on frequency through the dielectric constant. 
This parameter is often associated with the field within 
an empty cavity in a liquid [2 if. An expression recently 
derived by us for this property [30j was used in the calcu- 
lations: f d (u) = [7(e(w) + l) 2 +8e(w)]/[12e(w)(2e(w) + l)]. 

Since the polarizability of many organic substances is 
close to ao, e = <Tq/16, the parameter of dipolar density 
of the solutes [Eq. (fl6)) ] was taken in the form 



2/0 M 



+ ^kK) 2 (i-«W) 



m, (A2) 



where (m-o) 2 = /3(M 2 )/oq is a reduced effective dipole, 
M is the entire dipole moment of the protein- water clus- 
ter. 

Simulations of hydrated plastocyanin were reported 
previously [l8[. The presently used data represent 
the same simulation protocol applied to the oxidized (to- 
tal charge —8) state of plastocyanin extended to a larger 
number of waters in the simulation box, N w = 21076. For 
plastocyanin calculations M represents the total dipole 
of the protein and water shell extending 20 A from 
the protein surface into the bulk. This latter magni- 
tude was added to the effective radius of the protein 
listed in Table Q] to obtain the effective radius of the 
water/protein cluster. The averaged square of the pro- 
tein/water dipole calculated from the simulation trajec- 
tory was (M 2 ) = 1.44 x 10 6 D 2 . The response function 
was obtained as a Laplace-Fourier transform of the 



*(*) = !> 



-th 



(A3) 



where A, = {0.84,0.11,0.05} and n = {0.14,1790,6.3} 
ps. 

Other solute parameters used in the calculations are 
listed in Table HI the hard sphere diameter of water was 
taken at the value of a = 2.87 A, and the inertial param- 
eter p' in Eq. (J2UJ) was set at the value of v' = 0.1 [27f. 
The rotational relaxation times of the solutes were taken 
at To = 50 ps for trehalose and To ~ 3 ns for the two pro- 
teins. The former number is consistent with the second 
relaxation process extracted from the dielectric response 
and simulations of hydrated saccharides [50(, while the 
latter is typical for rotational dynamics of proteins [4l| . 

The calculation of the solute dipole component of the 
dielectric response simplifies in the continuum limit when 
the integral in Eq. (fTg)) loses the dependence on frequency 
and reduces to Eq. (f2"Tj) . This integral depends on two 
parameters, the volume fraction 770 and the reduced ge- 
ometry parameter r — 1/2 + cr/(2<7o), when the hard- 
spheres approximation is used for the density structure 
factor So(k,r]o). The range 0.5 < r < 1 covers most 
problems of interest. Numerical integration of Eq. (f2Tj) 
with the PY density structure factor [25[ was done in 
this range of r- values and volume fractions in the range 
< ?7o < 0.3. The numerical results were interpolated 
with the polynomial function 

I( m , r) = 0(770) + b{ m )r 2 + c(77o)r 4 + d(rio)r 6 , (A4) 



where 



1 + 0.225% 



aim) 

6(770) = -9.694770 - 18.57277 2 
c(?7o) = 6.987t7o + 38.913T7 2 - 
d(?]o) = 2.IO8770 - 16.5707?^ - 



7.72677^ - 13.8057?;] 
f 16.642770* 
5.940^ 
10.00777^ 



(A5) 
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The expansion in even powers in r in Eq. (|A4[) is dictated the density expansion of the polynomial coefficients has 
by the symmetry of the density structure factor [25[, and been chosen to justify the ideal-solution limit 1(0, r) = 1. 
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